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This paper provides central limit theorems for the wavelet packet decomposition of stationary band- 
limited random processes. The asymptotic analysis is performed for the sequences of the wavelet packet 
coefficients returned at the nodes of any given path of the M -band wavelet packet decomposition tree. It is 
shown that if the input process is centred and strictly stationary, these sequences converge in distribution 
to white Gaussian processes when the resolution level increases, provided that the decomposition filters 
satisfy a suitable property of regularity. For any given path, the variance of the limit white Gaussian 
0^ . process directly relates to the value of the input process power spectral density at a specific frequency. 
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5h ' I. Introduction 

This paper addresses the statistical properties of the A/-Band Discrete Wavelet Packet Transform, 
hereafter abbreviated as M-DWPT. In [1], asymptotic analysis is given for the correlation structure and 
the distribution of the M-Band wavelet packet coefficients of stationary random processes. The limit 
autocorrelation functions and distributions are shown to be the same for every M-DWPT path. This 
seems to be a paradox because the M-DWPT paths are characterised by several sequences of wavelet 
filters. Two arbitrary sequences are different, and thus, do not have the same properties. In addition, the 
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results presented in [1] seem to be in contradiction with those stated in [2] concerning the autocorrelation 
functions of the standard discrete wavelet packet transform (M-DWPT with M = 2) of wide-sense 
stationary random processes. The results presented in [2] highlight that the limit variance of the wavelet 
packet coefficients does depend on the path followed in the wavelet packet tree. 

In fact, as shown below, the limit results in [1, Corollary 5, Proposition 12] apply to only one path of 
the M-DWPT, namely the standard approximation path (only low-pass filters are used in this path). The 
same holds true for [3, Proposition 7], which extends [1, Corollary 5] to the case of the dual-tree M- 
DWPT Actually, the limit analysis of the autocorrelation and distributions of the M-DWPT coefficients 
is more intricate than presented in [1] and [3] because, as shown below, this analysis depends on the 
path chosen and the wavelet filters used for decomposing the input random process. This analysis is 
presented for the Shannon M-DWPT filters and standard families of paraunitary filters that converge to 
the Shannon paraunitary filters. 

II. Preliminary results 

A. General formulas on the M-DWPT 

In what follows, j and M are natural numbers and M ^ 2. An M-DWPT is performed by using 
wavelet paraunitary filters with impulse responses h m , m = 0, 1, 2, . . . , M — 1. For further details about 
the computation and the properties of M-DWPT filters, the reader is asked to refer to [4]. Let 

H m (u) = —= ^2 h m[i] exp (-iluj) , (1) 

and $ be a function such that {t^ : k G Z} is an orthonormal system of L 2 (R), where : t i — > 
<3?(i — k). Let U be the closure of the space spanned by this orthonormal system. 

The M-DWPT decomposition of the function space U involves splitting U into M orthogonal sub- 
spaces (an easy extension [5, Lemma 10.5.1] established for the standard DWPT): 

M-l 

U = W i,™> (2) 

m=0 

and recursively applying the following splitting 

M-l 

W i.™ = ® W i+l,M«+m, (3) 
m=0 

for every natural number j and every n = 0, 1, 2, . . . , M J — 1. In this decomposition, the wavelet packet 
space W j n is the closure of the space spanned by the orthonormal set of the wavelet packet functions 
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{Wj,n,k k <E Z} whose Fourier transforms are given by 

^Wj,n,fc(w) = exp (-iMiku) FWj, n (u), (4) 

with 

TWj, n {u) = M3' 2 FW n {Miuo), (5) 
where the sequence W n , n = 0, 1, 2, . . . , is recursively defined by 

FW m (oj) = H m (u/M)F<S>(u/M), (6) 

and 

FW Mn +m{w) = H m {u/M)TW n {u/M), (7) 

for m = 0, 1, 2, . . . , M - 1 and n € N. 

Note that the function <I> in Eq. (6) is not necessarily the standard scaling function associated with 
the low-pass filter h (see [5, Lemma 10.5.1] for more details). If <I> is this scaling function, we have 
W = $ in Eq. (6). 

B. M-ary representations of the M-DWPT paths 

A given wavelet packet path V is described by a sequence of nested functional subspaces: V = 
(U, {Wj jT j(j)}j 6 pj)> where Wj j7l (j) C Wj_i |tl ^_i), with n(0) = by convention and 

n(j) =Mn(i-l)+TO j; (8) 

for j ^ 1, where G {0, 1, . . . , M — 1}. Therefore, the shift parameter is 

j 

n(j) = ™ e MJ- e G {0, 1, ... , M J — 1} (9) 
e=i 

at every resolution level j. By construction, each ~Wj >n (j) is obtained by recursively decomposing U via 
a particular sequence of filters (h mt )e=i,2,...,j where each rri£ belongs to {0, 1, . . . , M — 1}. Thus, path V 
can be associated with a unique M-ary sequence A = (m^g^ of elements of {0, 1, . . . , M — 1}. From 
now on, any given M-DWPT decomposition path will be represented by an M-ary sequence A. Since 
the shift parameter n depends on j and A via Eq. (9), the notation n = n\(j) will hereafter be used 
to indicate this dependence if required. Therefore, an M-ary sequence A associated with an M-DWPT 
decomposition path specifies a unique sequence {^j t n x (j))jeN of wavelet packets. Now we have: 

Lemma 1: Let V = (U, {^j,nx(j)}j&) rje some P at h of the M-DWPT decomposition tree. If the 
shift parameter n\ is a bounded function of j, then A is the null sequence Ao = (0, 0, . . .). 
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Proof: If there exists some jo sucn that n\(jo) / 0, then it follows from Eq. (9) that n\(j) ^ M J j0 
for every j ^ jo- Therefore, nA cannot be upper-bounded by a constant. ■ 

Remark 1: Lemma 1 states that the approximation path (associated with the null sequence) is the unique 
M-DWPT decomposition path for which the shift parameter n is bounded by a constant independent of 
j. As a consequence, the results established in [1, Corollary 5, Proposition 12] and [3, Proposition 7] 
concern the approximation path only, because the assumption that n is constant is made in this reference 
to derive the asymptotic analysis. 

The results given in the present paper depends on the path chosen in the M-DWPT tree. In fact, the 
M-ary representation of the M-DWPT paths plays an important role throughout because it allows a 
complete path characterization. The following lemma will prove useful in the sequel. 

Lemma 2: For n = n(j) given by Eq. (9), we have 



i[H mt (- 



I=i 



F*(inr)- ( 10 ) 



Proof: An easy extension of [2, Lemma 1]. 



C. Shannon M-DWPT and the Paley -Wiener space of it band-limited functions 

The Shannon M-DWPT filters are hereafter denoted hf n for m = 0, 1, . . . , M — 1. These filters are 
ideal low-pass, band-pass and high-pass filters. We have 

fl* M = £ l Af > - 27rf), (11) 
eel, 

where Ik denotes the indicator function of a given set K: 1k(x) = 1 if x G K and \k{x) = 0, 



otherwise, and A r 



(TO+l)7r TO7T 

M ' M 



u 



m-K (m+l)n 
M ' M 



. The scaling function <& associated with these 



filters is defined for every t G R by <I> s (t) = sinc(£) = sin(irt)/irt with <1> S (0) = 1. The Fourier transform 
of this scaling function is 

^ s = l[_ w , ff] . (12) 

The closure U s of the space spanned by the orthonormal system {rfc$ s : k € Z} is then the Paley- 
Wiener (PW) space of those elements of L 2 (IR) that are it band-limited in the sense that their Fourier 
transform is supported within [— 7r, tt]. 

Let X be any band-limited Wide-Sense Stationary (WSS) random process whose spectrum is supported 
within [— 7r, tt]. We have (see [2, Appendix D]) 

X[k]= [ X(t)$ s (t- k)dt, (13) 
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so that U s is the natural representation space of such a process. Any M-DWPT of X can thus be 
initialized with the samples X[k], k £ Z. 

Now, let us consider the Shannon M-DWPT of the PW space U s . The wavelet packet functions Wj n 
of this M-DWPT can be computed by means of Eqs. (5), (6) and (7), by setting $ = <1> S and H m = H^, 
m = 0, 1,...,M — 1. The Fourier transforms of these wavelet packet functions are given by proposition 
1 below, which extends [6, Proposition 8.2, p. 328] since the latter follows from the former with M = 2. 



Proposition 1: For every non-negative integer j and every n € {0, . . . , M J — 1}, 

FW 5 = M J / 2 1a 

J yy J,n ly± 'A^g^j! 

where, for any non-negative integer k, 



(14) 



{k + l)ir kn 
W ' ~~M~3 



U 



kn {k + 1)tt 



(15) 



and G is defined by recursively setting, for m = 0, 1, . . . , M — 1 and I = 0, 1, 2, 

MG(£)+m if G(£)is even, 

MG(£) -m + M-l if G{£) is odd. 



G(0) = and 



G{Ml + m) 



(16) 



In the rest of the paper, we set, for any pair (j, k) of non-negative integers, 

kir (k + l)ir 



Mi ' Mi 



(17) 



III. Asymptotic analysis for the autocorrelation functions of the M-DWPT of 

SECOND-ORDER WSS RANDOM PROCESSES 

Let X denote a centred second-order real random process assumed to be continuous in quadratic mean. 
The autocorrelation function of X, denoted by R, is defined by R(t, s) = E[X(t)X(s)]. The projection 
of X on Wj „ yields a sequence of random variables, the wavelet packet coefficients of X: 



j, n [k] = [ X(t)W jtn , k (t)dt, k£Z, 
JR 



provided that the integral 



.// 



R(t, s)W j ^ k (t)W j ^ k (s)dtds 



(18) 



(19) 



exists, which will be assumed in the rest of the paper since commonly used wavelet functions are 
compactly supported or have sufficiently fast decay. The sequence given by Eq. (18) defines the discrete 
random process c Jiri = (cj^ n \k})k£i of the wavelet packet coefficients of X at any resolution level j and 
for any shift parameter n € {0, 1, ... , M J — 1}. 
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A. Problem formulation 

Let Rj :Tl stand for the autocorrelation function of the random process Cj :Tl . We have 

R jtn [k,e\=E[cj, n [k]cj, n [£\]= [[ R(t,8)W jtnt k(t)W jtn 4s)dtds. (20) 

J JR 2 

If X is WSS, we write R(t, s) = R(t — s) with some usual and slight abuse of language. From Eq. 
(20), it follows that 

Rj, n \ k A = II R(t)W jinjk (t + s)W jin/ {s)dtds. (21) 

J JR 2 

In the sequel, the spectrum 7 of X, that is, the Fourier transform of R, is assumed to exist. By using 
Fubini's theorem and Parseval's equality, we can proceed as in [2, Appendix C] to derive from Eqs. (4), 
(5) and (21) that c^ n is WSS. For any k,£ G Z, and with the same abuse of language as above, the value 
Rj,n[k,£\ of the autocorrelation function of the discrete random process Cj j7l is Rj, n [k — £] with 

R] ' n[k] = ^ ^(^WHI'exp^M^- (22) 
Let us assume that 7 G L°°(R) and is continuous at 0. These two assumptions have two easy 
consequences. First, the integrand on the right hand side (rhs) of Eq. (22) is integrable since its absolute 
value is upper-bounded by ||7||oo|^ 7 W / n(-)| 2 ' whose integral equals ||'y ||ooi second, the limit of 7(7^7) 
is 7(0) when j tends to 00. Therefore, for every given natural number n, it follows from Lebesgue's 
dominated convergence theorem applied to Eq. (22) that 

lim Rj, n [k] = 7T [ 7(0)|JP"W„(u;)| 2 exp (ikw) du , 

= 7(0) / W n (t)W n (t - k)dt = 7 (0)<J[A;], (23) 

JR 

1 if m = 0, 



if m^O. 



where 5 is the standard Kronecker symbol defined for k G Z by 5[k] = 

The result thus obtained is that given in [1, Corollary 5]. 

From Lemma 1, we distinguish two cases, for any given M-DWPT path V = (U, {Wj nx Q\}j e jqj. 
First, if n\ is a constant function of j, then A is the null sequence Ao, and thus, V is the approximation 
path. In this case, the shift parameter n\(j) is at each resolution level j and the M-DWPT of X 
through path V = V\ Q consists of an infinite sequence of low-pass filters. The decorrelation is then 
guaranteed by Eq. (23) (see also [1, Corollary 5]). The second case is that of a function n\ which cannot 
be upper-bounded by a constant independent of j when j tends to infinity 1 . In such cases where n\ is 

Example: for the sequence A = (1, 1, . . .), we have n\{j) = M 3 — 1 . The nodes (j, M J — 1) are those of the path located 
at the extreme rhs of the M-DWPT decomposition tree. 
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not a constant function of j, the asymptotic decorrelation of the M-DWPT coefficients at node (j, n\(j)) 
when j tends to oo is no longer a mere consequence of Eq. (23) since Lebesgue's dominated convergence 
theorem does not apply. To proceed, we then write Eq. (22) in the form 

RjAk] = [ lH\FW j>n (uj)\ 2 exp UM'kw) dco. (24) 

This equality derives from Eq. (22) after a change of variable and by taking into account Eq. (5). 

The purpose of the next section is then to analyse the behaviour of Rj iU in the case of the Shannon 
filters and some families of filters that converge to the Shannon filters. 

B. Asymptotic decorrelation achieved by the Shannon M-DWPT 

Let A = (mejem be an M-ary sequence of elements of {0, 1, . . . , M — 1}. Consider the Shannon M- 
DWPT, that is, the decomposition of U s associated with the Shannon M-DWPT filters (/^) m =o,i,...,M-i- 
Let V\ = (U s , {W| n A (j)}jeN) be the path associated with A in the Shannon M-DWPT decomposition 
tree. It follows from proposition 1 that the support of is Aj jPx (j^, where p\(j) = G(n\(j)). For 

j G N, the sets are nested closed intervals whose diameters tend to 0. Therefore, their intersection 

contains only one point oj\. It then follows from (15) that 

J^ + OO M-? 

Let X be some centred second-order WSS random process, continuous in quadratic mean, with 
spectrum 7. The autocorrelation function Rj n resulting from the projection of X on W| n derives from 
Eq. (24) and is given by 

Rj,n[k] = ^ / 7 (^)l^ n ( w )| 2 exp (iMlku) dw. (26) 

From Eqs. (14) and (26) and by taking into account that 7 is even, as the Fourier transform of the even 
function R, it follows that 

R$ n [k] = / 7(0;) cos (M^kco)dw. (27) 

vr Ja+ 

JyP 

where A^ p is given by Eq. (17) and p = G{n). When X satisfies some additional assumptions, the 
following theorem 1 states that the Shannon M-DWPT of X yields coefficients that tend to be decorrelated 
when j tends to infinity. One of these additional assumptions is that X is band-limited in the sense that 
its spectrum is supported within [— 7r,7r]. When M = 2, theorem 1 is equivalent to [2, Proposition 1]. 

Theorem 1: Let X be a centred second-order WSS random process, continuous in quadratic mean. 
Assume that the spectrum 7 of X is an element of L°°(R) and is supported within [— 7r,7r]. Let A = 
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(mg)eeN be an M-ary sequence of elements of {0, 1, . . . , M — 1} and V\ = (U s , {W| n ^ ;)}jeN) be the 
Shannon M-DWPT decomposition path associated with A. 
If the spectrum 7 of X is continuous at point uj\, then 




uniformly in k G Z, where i? s , s is the autocorrelation function of the coefficients resulting from the 



Proof: The proof is an easy generalisation of that of [2, Proposition 1], which concerns the standard 
wavelet packet transform (M = 2). The key point of the proof is proposition 1 above, which makes it 



The foregoing theorem is mainly of theoretical interest since the Shannon M-DWPT filters have infinite 
supports and are not really suitable for practical purpose. In order to obtain a result of the same type 
for filters of practical interest, the M-DWPT is now assumed to be performed by using decomposition 
filters of order r, h\n, m = 0, 1, . . . , M — 1, such that 



The parameter r is called the order of the M-DWPT filters. According to [7], the Daubechies filters 
satisfy Eq. (29) for M = 2 when r is the number of vanishing moments of the Daubechies wavelet 
function; according to [8], Battle-Lemarie filters also satisfy Eq. (29) for M = 2 when r is the spline 
order of the Battle-Lemarie scaling function. The existence of such families for M > 2 remains an open 
issue to address in forthcoming work. However, it seems reasonable to expect that general M-DWPT 
filters of the Daubechies or Battle-Lemarie type converge to the Shannon filters in the sense given above. 

Theorem 2: Let X be a centred second-order WSS random process, continuous in quadratic mean. 
Assume that the spectrum 7 of X is an element of L°°(R) and is supported within [— tt, it]. Assume that 
the M-DWPT of the PW space U s is achieved by using decomposition filters hm, m = 0,1, . . . , M — 1, 
satisfying Eq. (29). 

For every natural number j and every n = 0, 1, . . . ,M J — 1, let R - n stand for the autocorrelation 

[rl 

function of the wavelet packet coefficients of X with respect to the packet 1 . We have 



projection of X on MV^ 



possible to compute Eq. (25). 



limtfH=tfS (a.e.). 



(29) 



lim R^[k]=Rl n [k], 



(30) 



uniformly in k G Z and n, where 



is given by Eq. (27). 
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Remark 2: Since the Af-DWPT concerns the space U s , we have 



3 



LO 



3 



W 

i=i 



LO 

'"' y M3 +l - 



(31) 



(32) 



where n is given by Eq. (9). These equations straightforwardly derive from Eq. (10) of lemma 2. From 
Eqs. (5), (29), (31) and (32), we obtain, for every given natural number j, that 



lim rwfl = TWl n (a.e.), 



r— >+oo • ?,ri " ' h 

uniformly in n. The three equalities above will prove useful below. 



(33) 



\r\ 

Proof: (of theorem 2). The autocorrelation function Ry n is given by Eq. (24) and is equal to 
In addition, we have 

\R^[k] ~ *?,„[*] I < ^ jf I7HI | l-™$MI 2 " dLJ, (35) 

where Rj n is given by Eq. (27). From Eqs. (5), (31) and (32), and by taking into acount that 

and \Hf n (lo)\ are less than or equal to 1 (due to the paraunitarity of the M-DWPT filters), we obtain 



FWJI{lo)\ 2 -\FWUuo) 



< 2M 3 



T&{lo) 



(36) 



The results derives from Eqs. (33), (35), (36) and Lebesgue's dominated convergence theorem. 



IV. Central limit theorems 
We now consider a real random process X that has finite cumulants and polyspectra. Denote by 

cum(i, si, 32, ■ ■ ■ , s N ) = cum{X(t),X(si),X(s 2 ), X(s N )}, (37) 

the cumulant of order N + 1 of X. The above cumulant is hereafter assumed to belong to L 2 (M. N+1 ) 
and to be finite for any natural number N (see [9, Proposition 1] for a discussion about the existence of 
this cumulant). The cumulant of order N + 1 of the random process Cj n has the integral form given by 
(see [9, Proposition 1]): 

cum j)n [A;,^i,...,^v] = cum { c j,n[ k ] c 3,n[^} ■■■ c j,n[^N]} 

dtdsi . . .ds N cum(t,s 1 ,s 2 , ■ ■ . , SAr)Wj>,fe(*)^>A( s i) ■ ■ ■ w j,n,e N (s N )- (38) 



/ 
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Assume that X is strictly stationary so that cum(i, t + t\, t + t2, ■ ■ ■ , t + tjy) = cum(ti, t2, ■ ■ ■ , tjy), 
then Cj :n is a strictly stationary random process with cumulants cum j in [k, k + k±, k + fe, . . . , k + Ivn] = 
cumj n [ki, k2, • • • , fejy]. Assume also that X has a polyspectrum 7^(0^, 0J2, • • • , c^v) £ L°°(M. N ) for 
every natural number N and every (ui,U2, ■ ■ ■ ,wn) £ ^ N ■ The polyspectrum is the Fourier transform 
of the cumulant cum(ti, t2, ■ ■ ■ ,tN). When N = 1, 71 is the spectrum of X and is simply denoted 7 as 
in section III. Then after some routine algebra, Eq. (38) reduces: 

M -j(N-l)/2 p 

cum jn [ki,...,k N ] = — — / dun ■ ■ • dwArexp (-i{kiw\ + . . . + k N u; N )) 

(27T) JV JR" 

• • • 7jv (-u;iM-J, . . . , -u N M-j)^W n (- Ul - ... - uo N )TW n (^) . . . TW n {u N ). (39) 

If n is a bounded function of j, it follows from Lebesgue's dominated convergence theorem that, for 
any natural number N > 1, curr\j^ n [ki, k2, ■ ■ ■ , fcjv] tends to uniformly in fci, fo, . . . , fcjv when j tends 
to 00. This is a consequence of [1, Proposition 11]. On the other hand, if n cannot be upper-bounded 
by a constant independent with j, the situation is similar to that discussed in section III- A: the shift 
parameter n depends on j and Lebesgue's dominated convergence theorem does not apply to Eq. (39) to 
prove the vanishing behaviour of the cumulants. The vanishing behaviour of the cumulants are hereafter 
given by theorems 3 and 4. 

By taking into account Eq. (5), Eq. (39) can also be written, with an easy change of variables: 

cumj jri [fci, ...,k N ]= 1 N j dwi . ..div N exp {-iM j (kiUJi + . . . + k N iv N )) 

7jv(-^i, . . . , -u N )FWj >n (-ui - ... - uj N )TW jtn (uJi) . . . FW j>n (u N ). (40) 

Theorem 3: Let X be a centred second-order strictly stationary random process, continuous in quadratic 
mean. Assume that the polyspectrum 7^ of X is an element of L°°(R N ) for any N ^ 1 and that the 
spectrum 7 is supported within [— ir, tt\. If N > 1, then we have, uniformly in fci, k2, ■ ■ ■ , fc/v, 

lim cumf [A;i, fo, . . . , fcjv] = 0. (41) 

Proof: When the wavelet packet functions are the functions Wj n , it follows from Eqs. (14) and (40) 
that the cumulant cum| [fci, k2, ■ ■ ■ , /cat] of the discrete random process returned at node (j, n) by the 
Shannon M-DWPT of X satisfies 

|cumf n [fci, k 2 , ■ ■ ■ ,kw]\ ^ jz , N — — / dwidw 2 ...dujN (42) 

(27T) JV J A n 

J,P 

where A^ p = Aj iP x Aj iP x . . . x Aj tP and p = G(n). 
N times 
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According to Eq. (15), J A dco = 2t:/MK Therefore, we obtain 

|cumf jn [fci, k 2 , ...,k N ]\< bivllooM-^- 1 )/ 2 . (43) 

Given any natural number N > 1, the rhs of the latter inequality does not depend on n, k±,. . . ,k>N and 
vanishes when j tends to oo, which completes the proof. ■ 

Corollary 1: With the same assumptions and notation as those of theorems 1 and 3, assume that 7 is 
continuous at lu\. Then, when j tend to infinity, the sequence ( c j„ A (j)) converges in the following 
distributional sense to a white Gaussian process with variance 7(0;^): for every x 6 R N and every e > 0, 
there exists jo = jo(x,e) > such that, for every j ^ jo, the absolute value of the difference between 
the value at x of the probability distribution of the random vector 

and the value at x of the centred A-variate normal distribution J\f(0, j(u}\)In) with covariance matrix 
~/(uj\)In is less than e. 

Proof: A straightforward consequence of theorems 1 and 3. ■ 
Consider filters satisfying Eq. (29). Let A be an M-ary sequence of elements of {0, 1, ... ,M — 1}. 



\r] 

The following results describe the asymptotic distribution of the discrete random process c . (i) returned 



at node (j,n\(j)) when the resolution level j and the order r of the filters increase. 

Theorem 4: Let X be a centred second-order strictly stationary random process, continuous in quadratic 
mean. Assume that the poly spectrum 7^ of X is an element of L°°(R N ) for every natural number N ^ 1 
and that the spectrum 7 is supported within [— 7r,7r]. 

For every given natural number j and every n € {0, 1, ... , M- 7 — 1}, let cumj^ stand for the cumulant 

M 

of order N + 1 of the wavelet packet coefficients of X with respect to the packet ' n . 
We have, uniformly in n, k±, k 2 , • • • , 



lim cum [ - ] [k 1 ,k 2 ,...,k N ] = cum^ n [k 1 ,k 2 ,. . . ,k N ]. (44) 

Proof: By applying Eq. (40) to W^'\ and W^ n , we obtain 
cum^jjfci, . . . ,k N ] -cum^ n [k 1 ,...,k N ]\ < 1 „ | |-yjv| |oo / dwi ■ ■ ■ duj N 

' (27T) JV J r n 

rwflt-u! ... - u^rwfl^) . . . rwf)sm) - ^n(-wi • • • - u^rwj^) . . . rw? n (u N ) 

(45) 
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The integrand on the rhs of the second inequality above can now be upper-bounded by 

2M^ n+1 ^H s (lo 1 )^ s (uj 2 ) . . . $ s (u N ) (46) 

where we use Eqs. (5), (31), (32), and take into acount that |ilj^(a;)| and \Hf rli {uj)\ are less than or 
equal to 1. The upper-bound given by Eq. (46) is independent of r and integrable; its integral equals 
2M^ N+1 ^ 2 (2ir) N . By taking Eq. (33) into account, we derive from Lebegue's dominated convergence 
theorem that the upper bound in Eq. (45) tends to when r tends to +oo. ■ 

Corollary 2: With the same assumptions and notations as those of theorems 2 and 4, assume that 7 
is continuous at u\. Then, when j and r tend to infinity, the sequence (c^ ] converges in the 
following distributional sense to a white Gaussian process with variance j(oj\): for every x 6 M. N and 
every e > 0, there exists jo = jo(x, e) > and there exists tq = ro(x, jo, e) such that, for every j ^ jo 
and every r ^ tq, the absolute value of the difference between the value at x of the probability distribution 
of the random vector 

(c [r] r M,c [r] ( Jk 2 ],...,c [r] ( Jk N ]) 

and the value at x of the centred ./V-variate normal distribution A^(0, j(u)\)In) with covariance matrix 
~/(uj\)In is less than e. 

Proof: The result follows from Eqs. (28), (30), (41) and (44). ■ 
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